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We introduce a new scheme to analyze the kinetics of homogeneous nucleation in terms of a 
global order parameter. Our approach is based on the application of the internal degrees of freedom 
formalismcl to derive a kinetic equation of the Kramers type formulated for a global reaction coor- 
dinate. We provide explicit expressions for the quantities and coefficients involved in the process, 
suitable for simulation. In addition, our picture recovers in the quasi-stationary case the transition 
rate obtained from the method of reactive flux. The equation we present may provide a link be- 
tween theoretical approaches to homogeneous nucleation (generally formulated in terms of a kinetic 
equation of the Fokker-Planck type) and simulations (which mostly employ linear response theory). 
In this context, our scheme pttruides a theoretical framework to interpret and extend the results 
obtained in recent simulationsB □ 

I. INTRODUCTION 

Nucleation is the fundamental mechanism, of many phase transformations, and has been the subject of extensive 
theoretical and experimental investigations^ In the last years, computer simulations have also become a useful tool 
to study the kinetics of crystal nucleation. However, in order to analyze the results of simulations the most common 
theoretical treatments for nucleation processes are not generally employed. 

The kinetics of homogeneous nucleation is habitually described from a theoretical point of view by means of kinetic 
equations of the Fokker-Planck type focused on the size distribution of the clusters present in the system. These 
equations are formally analogous to the ones arising from Kramers reaction-rate theory,trQ since in fact nucleation is 
an activated process. The several approaches proposed to evaluate nucleation rates differ basically in the expressions 
for the drift and diffusion coefficients of the Fokker-Planck equation. The values of these coefficients are sometimes 
postulated, assumed as unspecified parameters, or left in terms of the addition and loss rates of one monomer to a 
cluster. To estimate these rates a proper microscopic kinetic model is required. These theories are customarily used 
to evaluate steady-state nucleation rates, and most of them assume the constancy of the diffusion coefficient. r-, 
Contrarily, two different ways to proceed exist for simulations. On the one hand, there is a type of simulation,El 
which can be called " straightforward" , whose method consist of simply supercooling the system and then waiting for 
the appearance of nuclei. Since the nucleation barrier decreases as 1/(AT) 2 , where AT is the degree of supercooling, 
this type of simulation based on observing the spontaneous formation of nuclei requires experimentally inaccessible 
supercoolings. For more moderate and realistic conditions the nucleation barrier is so large that a spontaneous crossing 
is extremely unlikely during the typical times of a simulation. 

Recently, a different type of simulation, allowing the study of nucleation under more realistic undercoolings, has 
been introduced by ten Wolde, Ruiz-Montero and Frenkel in refs. |[^. Their method is based on the separation of. 
the simulation into two parts: the evaluation of the nucleation barrier following the Van Duijneveldt^Frenkel schemed 
and the subsequent calculation of the nucleation rate by means of the Bennett-Chandler pictureJi3~E3 based on linear 
response theory. The macroscopic state of the simulated system is described by means of a reaction coordinate 
connecting the initial and final phases and sensitive to the global degree of crystallinity Loathe system. Concretely, 
Frenkel and co-workers use a set of bond-order parameters introduced by Steinhardt et al&3 

One may wonder why they use linear response theory instead of appealing to the usual nucleating theories or even 
$_i ' to the prototypical theory of activated processes: the Kramers theory. In fact, in their simulations&tj they report on a 
diffusive behavior that seems to corroborate the validity of Kramers picture. Additionally, they also employ Kramers 
equation in order to obtain reasonable estimates of the time evolution of the system. 

In reference the authors themselves give the reasons why they do not use Kramers theory. First, they assert that 
it was originally formulated for the single particle diffusion problem and therefore its validity is not guaranteed in the 
case when the reaction coordinate is a global one. On the other hand, they point out that the diffusion coefficient, 
assumed to be constant in Kramers picture, may depend on the value of the reaction coordinate, and moreover its 
expression may be non-trivial. Finally, they affirm that a Fokker-Planck description of nucleation processes does not 
take into account hydrodynamic effects that may be very important in this type of barrier-crossing processes. 



. — i 

X 



1 



Our purpose in this paper is to introduce a theoretical description of homogeneous nucleation in terms of a global 
order parameter like the one used by Frenkel and co-wockers in their simulations. The treatment we present is based 
on the kinetic equation introduced in a previous paperta and derived from the application of the internal degrees of 
freedom formalism. EH13 We want to show how the scheme we propose provides a complete description of the process, 
it includes the transition rate expression of the method of reactive flux as a quasi-stationary limit, and it overcomes 
the reservations exposed by Frenkel et al. 

The paper is structured as follows. In section II, we outline and briefly review the method of reactive flux from 
which expressions used by Frenkel's group in their simulations are derived. In section III, we introduce and extend the 
internal degrees of freedom formalism to the situation where the internal coordinate is a non-local order parameter. 
Within this framework, we derive a Kramers-like equation, which in section IV is particularized to the quasi-stationary 
case to prove that it recovers the results of the method of reactive flux. In section V, we describe the way to proceed 
beyond the quasi-stationary limit. Finally, in the last section we summarize and comment on the main results of the 
theory we propose. 



II. METHOD OF REACTIVE FLUX. 



Our aim in this section is to present a brief review of the derivation, the limitations, and the scope of applicability 
of the method of reactive flux to evaluate reaction rates. To this end, let us consider a macroscopic system whose 
state is characterized by a reaction coordinate Q. This system can be in two different states, A and B, separated 
by an energy barrier whose maximum is located at position Qq . The state of the system, i.e. the side of the barrier 
where it is located, can be characterized by the functions 

n A (t) = 9[Q - Q(t)] (1) 



n B [t) = 6[Q{t) - Qo] (2) 

where 9 is the step function, and n A (t), n B (t) obviously satisfy the condition 

n A {t) + n B {t) = 1 (3) 

The probability that our system is in state A at time t will be given by the non-equilibrium average of the characteristic 
function n A (t) 

P A (t) = (n A (t)) (4) 

The starting point of this theory is to postulate the validity of the following phenomenological equation describing 
the dynamics of the populations in states A and B 

- k BA P B {t) - k AB P A (t) (5) 

at 

where P B (t) = 1 — P A (t), and k AB and k BA are the forward and backward rates, respectively. These rates are 
assumed to be constant and satisfy the detailed balance condition 



k AB PJ_ 

r>eq 
A 



- peg ( 6 ) 
k BA Pa 



In the spirit of Onsager's regression hypothesis, we assume that the relaxation of AP A (t) = P A (t) — P^ 1 from an 
initial non-equilibrium deviation AP A (0) follows the same exponential decay as the equilibrium correlation function 
of the fluctuations 

AP4t) = {An A (0)An A (t)) eq = _ xt 

AP A (0) (An A (0)An A (0)) eq [ 1 

where An A (t) = n A (t) — (n A ) eq , the symbol (. . .) eq represents an equilibrium average, and the relaxation rate A is 
given by 

A = k AB + k BA (8) 
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From these equations one can obtain the-ibllowing formula for the transition rate, employed by ten Wolde, Ruiz- 
Montero and Frenkel in their simulationsoQ 



(h B (0)n B (t)) eq _ (Q(0)6[Q(0)-Q ]9[Q(t)-Q }) eq 



(nA)eq {nA)eq 



(9) 



This expression is not valid for very short times (because the relaxation of the system cannot be exponential at t = 0). 
Additionally, it is also restricted to times t satisfying the condition t << 1/A, since during the derivation it is assumed 
that e _A * ~ 1. Moreover, we notice that although k AB is a time independent rate, in expression (||) it is equated to a 
correlation function depending explicitly on time. Therefore, the latter equation is only correct for times long enough 
that the correlations have reached a plateau value. 
On the other hand, taking into account the identity 

n B {t) -n B (0) = j'dt'^jP- (10) 

and since the equilibrium average (h B n B (0)) eq vanishes, one can rewrite eq. (Q) as follows 

kAB = f (Mo)*^* 



{nA]eq 

As asserted before, this expression is only valid for times long enough that the correlations have achieved a plateau. 
For such values of the time, it is then meaningful to assume that the velocity autocorrelation function has decayed 
to zero, so that we can replace the upper limit of the previous integral by infinity. Therefore, the final expression for 
the transition rate yields 

k AB = / (n B {0)n B (t)) eq dt = — — / Q(0 )S[Q(0 ) - Q ]n B (t )} eq dt (12) 

{nA)eq Jo {nA)eq J 

This expression will be recovered by using the treatment we introduce in the next section. 



III. NON-EQUILIBRIUM THERMODYNAMICS APPROACH. INTERNAL DEGREES OF FREEDOM 

The key element of the approach we propose is the formulation of a kinetic equation to describe homogeneous 
nucieation in terms of the order parameter Q characterizing_the global degree of crystallinity of the system. This 
equation has already been introduced in a previous paperJU and we will devote this section to an extended and 
complete derivation of this equation. The procedure is based on the application of the internal degrees of freedom 
formalism, which constitutes a straightforward method to obtain kinetic equations for any adequate internal coordinate 
or degree of freedom of the system. In particular, we use this formalism considering Q as the system internal coordinate. 

Until now, in all of previous papers where this formalism has been successfully applied Jd'E^nEil the system was 
always described by means of a local internal coordinate identifying each system constituent independently. As a 
result, a Kramers-like equation was obtained. However, we are now interested in a kinetic description in terms of a 
single global order parameter Q, characterizing the whole system. Therefore, we need to extend the internal degrees 
of freedom formalism to allow a description in terms of a global internal coordinate. In fact, we face the same problem 
pointed out by Frenkel and co-workersEj about the validity and the necessity of an extension of Kramers picture when 
the reaction coordinate is a global one. .— . .— . 

From the point of view adopted in previous papersJijEa one worked with a system of N identical particles. The 
internal coordinate characterized the state of each particle, which evolved independently without interacting with the 
remaining particles. The global state of the system was obtained at each instant by means of a statistical recount of 
the individual states of the N particles. 

In our case, however, we are interested in a contracted description of the system based on a single internal coordinate 
that characterizes it in a global way. Proceeding by analogy, we consider N system replicas, identical in the sense that 
all of them are compatible with a given macroscopic state (identified for example by the total number of particles, the 
pressure, and the temperature). This set of replicas constitutes a non-equilibrium ensemble, defining the "system" we 
are going to deal with. Since we are working with an ensemble, we will be able to relate and obtain statistically the 
properties of a single replica from the behavior of the whole ensemble. Denoting the probability density by p(Q, i), the 
quantity p(Q, t)dQ represents the fraction of replicas with a value of the degree of crystallinity between (Q, Q + dQ), 
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at the instant t. For N large enough, the latter expression could be interpreted as the probability that a single replica 
of the ensemble has a value of the internal coordinate between (Q, Q + dQ). 

Therefore, the value of a quantity A evaluated for a single replica could be assimilated with the expected value of 
this observable in the non-equilibrium ensemble characterized by p(Q,t). That is 

A(t)\ repUca (A(t)) ensemUe = / Ap(Q,t)dQ (13) 

JQi 

Again, in order to derive our kinetic equation, our starting point will be the Gibbs equation for the entropy variation 
associated with a diffusion process through the energy barrier that separates the states Q\ (the value of the order 
parameter in the metastable phase) and Q 2 (corresponding to the new stable phase) 

SS = - i f ' p(Q, t)S P (Q, t)dQ (14) 
1 JQi 

Here 5S is the variation of entropy per replica, T is the temperature (taken as constant) and p(Q,t) is the conjugate 
chemical potential. The latter equation corresponds to the variation of entropy in a diffusion process, with the 
assumption that the internal energy and the volume of the system are constant. 

It is worth pointing out that the Gibbs equation constitutes a reformulation of the Gibbs entropy postulate. Indeed, 
from the fact that we are working with system replicas that obviously do not interact with each other, the expression 
for the chemical potential is that for an ideal system 

p{Q,t) = k B Tbxp{Q,t) + $(Q) (15) 

where $((3) is the energy barrier to be surpassed by a replica in order to achieve the new stable phase. Bearing in 
mind that at equilibrium the chemical potential becomes constant, the following expression for the energy barrier 
could be obtained from the latter equation 

*(Q) = p eq - k B T\n Peq {Q) (16) 

Therefore, we can rewrite the chemical potential as follows 

p(Q,t) = fx eq + fc B Tln^%^ (17) 



Peq{Q) 



Substitution of this equation in expression (|l4| ) yields 



SS = -k B 5p(Q,t)]n-?-dQ (18) 

J Peg 

which is just the expression for the variation of the entropy obtained from the Gibbs entropy postulate. Consequently, 
the equation we use as starting point in our derivation, namely eq. ([bi]), constitutes an alternative formulation of 
the Gibbs entropy postulate more proper for a non-equilibrium thermodynamics description. Therefore, it essentially 
represents the variation of 5* with respect to its local equilibrium value in terms of the distribution over the possible 
internal states associated with different values of the order parameter Q. 

The evolution of the number density p(Q,t) is governed by the continuity equation 

dp(Q,t) _ dJ(Q,t) 

8t dQ { ' 

where J(Q,t) is the diffusion current in the internal space. From eq. ( |l4| ) we can obtain the entropy production 
associated with this diffusion process 



(20) 



where we have employed eq. ( |19[ ) and a partial integration has been performed, assuming vanishing flux at initial and 
final states. 

The entropy production ( |20| ) has the usual form of a sum of flux-force pairs. Following the non-equilibrium ther- 
modynamics method, we may formulate the phcnomenological equation 
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m , t) . J&MjjS. (21 , 

where is the phcnomcnological cocfhcicnt, which may in general depend on the internal coordinate. To derive 

this expression, locality in the internal space has been assumed, for which only fluxes and forces with the same value 
of the internal coordinate are coupled. 

By substituting (^lj) in ( p^ ) and using eq. ( |l5| ) for the chemical potential, we obtain the following kinetic equation 
in the internal space 

% = £(mt& + KQ,t)g:p) (22) 



dt dQ V dQ dQ 

where 

- mm < 23 » 

and 

D(Q,t) = k B Tb{Q,t) (24) 

are the mobility and the diffusion coefficient in the internal space, respectively. 

Equation ( p2|) characterizes the evolution of the non-equilibrium physical ensemble, determining in a complete way 
the nucleation process in terms of the global order parameter Q. In order to fully specify this equation, we need an 
expression for the barrier $(Q) and for the phenomenological coefficient L(Q) relating the flux and its conjugate force 
in the internal space. 



The value of the energy barrier is given by eq. (16), and could be evaluated by means of simulations. In fact, this 
type of simulations for global order parameters involved in nucleation processes has already been performed by van 
Duijneveldt and Frenkel.EI 

On the other hand, to obtain an explicit expression for the phenomenological coefficient L(Q), we will adopt the 
fluctuating hydrodynamics framework.E3 From this point of view, one considers that the total current in the internal 
space is the result of two contributions: a systematic part J s and a random term J r including the fluctuations inherent 
in the diffusion process along the internal coordinate. One has 

J(Q,t) = J s (Q,t) + J r (Q,t) (25) 

We will assume that J r has zero mean and satisfies the fluctuation-dissipation theorem 

<J r (Q, t) J r (Q', t')) = 2k B L{Q)S{Q ~ Q')8{t - t') (26) 

Finally, from the latter equation one can derive the following Green-Kubo expression for the phenomenological coef- 
ficient 

1 pQ2 poo 

L (Q) = 1-1 d Q' dt(J r {Q',Q)J r (Q,t)) (27) 
k B Jq 1 Jo 

The objective of the two following sections will be to transform this expression into a formula suitable for computer 
simulations. First, we consider the quasi-stationary case to show how to recover explicitly the method of reactive flux 
expression commonly used to evaluate nucleation rates. Subsequently, we will analyze the general case. 



IV. THE TRANSITION RATE FOR THE QUASI-STATIONARY CASE. 



When the supercooling or the supersaturation is low enough, the energy barrier to surpass in the nucleation process 
is large compared with thermal energy. In this high barrier situation, the system achieves a quasi-stationary state 
characterized by a uniform current 

J(Q, t) = J(t) {9{Q - Qi) - 9{Q - Q 2 )} (28) 
and a chemical potential which equilibrates independently at each side of the barrier 
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M (Q, t) = p(Qi,t)9(Q - Q) + n(Q 2 , t)9(Q - Q ) , (29) 

where Qo represents the position of the top of the barrier and 9 is the step function. By substituting this equation in 
( |l5| ) one obtains the following expression for the probability density 

g(Q)-g(Qi) g(Q)-g(Q 2 ) 

p(Q,t) = p(Qi,t)e "bt e{Q a - Q) + P {Q 2 ,t)e ^^9(Q-Q ) (30) 
On the other hand, the diffusive current can be written in a more convenient way as 

J(Q,t) = -D(Q,t)e-^ k - T -^e^ k - T (31) 

By equating (|3l]) and (p8|), and by integrating over Q with the help of eq. (|29|), one arrives at the following expression 
for the nucleation rate 

J( t ) = ^(Qo) e -* o/kB T f e ,,/teT _ e ^/k B T\ (32) 

where the subindexes 0, 1 and 2 indicate that the corresponding function is evaluated at Qo, Q\ and Q2, respectively. 

Having obtained the nucleation rate expression for the quasi-stationary case, our goal now will be to show that it 
is equivalent to the method of reactive flux equation derived from linear response theory. An important difference 
between our formalism and the method followed in the simulations is that we are working with a non-equilibrium 
ensemble of replicas, while in simulations one makes reference to a single system or a single replica. The connecting 
point that links our scheme with the results referred to a single system is the probabilistic interpretation of the 
ensemble number density p(Q,t). 

Bearing this probabilistic interpretation in mind, we will see how to recover the transition rate expression (|l2|). The 
first step will be to transform the continuity equation into the phenomenological starting equation of the method of 
reactive flux. By integrating eq. ( [l9| ) from Q\ to Qo, and using the quasi-stationary condition (j2^) and the result 
( |32| ) one obtains 

4- f Q ° p(Q,t)dQ = -J(t) = ° {Q °l e -*°/*BT / m/kB T _ e „ l/kB T\ (33) 
at J Ql Q 2 - Qi V J 

where Jq" p(Q, t)dQ = Pa(£) IS j us t the probability that the system would be on the left hand side of the barrier, at 
time t. 

From ( p3| ) and (|24|) the diffusion coefficient expression at the top of the barrier is given by 

D(Qo,t) = k -^- (34) 
By using ( [l5| ) and (|30|), the density p{Qo,t) could alternatively be expressed as 

p(Qo,t) = e /*i/*BT e -» /kBT (35) 

On the other hand, from eq. (BOh and (FTs) it follows 



f Q* e -*/ kB T dQ e -*/ kB T dQ 

Here Pb represents the probability of being on the right hand side of the barrier. By introducing these two latter 
expressions into eq. (p3|), one arrives at 

= k BA (t)P B (t) - k AB {t)P A {t) (37) 

at 

where 

h (+\ L (Qo) k B 



G 



and 

, m L(Q )k B Q e-*"»TdQ 

kBA[t) ~ (Q 2 - Q,)P A {t) j£ e - W dQ (dJJ 

are the forward and backward rates, respectively. 

Therefore, we have just obtained the equation governing the evolution of the probability that the system would 
be on the metastable side of the barrier. At this point, we have to highlight two important differences with respect 
to the method of reactive flux. First, we have obtained this equation as a particular case (corresponding to a quasi- 
stationary situation) of a more general theory of wider applicability. Moreover, the rate coefficients in our expression 
are explicitly time-dependent, which does not happen in the formula of the method of reactive flux. Finally, it is 
important to remark that our rate coefficients fulfill the detailed balance condition 

eq 



k AB (t) _ pg 

A 



k BA {t) p e A q N(M 



as follows from eqs. (|T^), (|38| ) and (|3^). This condition guarantees that the flux vanishes at equilibrium. 

The last task to carry out is to compare the rate kAB{t) with the reactive flux formula (12). Within that framework, 
the expression (||) is identified with the rate Uab at times large enough so that the correlations have achieved a stable 
value characterized by a plateau. In our expression, for long times the probability Pa will approach and can be 
replaced by the equilibrium value P e A q — (riA)eq, an d therefore kAB will become time independent. 

On the other hand, the random component J r of the current is related to the fluctuations of P A with respect to 
the equilibrium value. Indeed, by integrating the continuity equation ( |l9| ) from Qi to Q one obtains 

,)S{QJ) = -J(Q,t) = -J s (Q,t) - J r (Q,t) (41) 



at 



where 



rQ 

N(Q,t) = / p(Q',t)dQ' (42) 



Notice that eq. ([id]) is valid for arbitrary times. In a very short time scale -microscopic time scaleci we can write 

^™ = -r(Q,t) (43) 

where the time derivative dSN ^Q-^ j s the instantaneous rate of change.^! Therefore 

(/(?.»)/(«,.»„ - ( a T ,0 "T'" »-« < 44 > 

with SN = N — N eq . In the quasi-stationary state, described by eq. (|2^), the quantity dSN gf is independent of the 
position Q' of the absorbing barrier. Hence, we can particularize for convenience the latter expression for Q' = Qq, 
and taking into account that 

^fA i (AlM(i)) . _| M m 

we can then rewrite the phenomenological coefficient (p7|) as 



k B L(Q ) = (Qi - Qi) / dt((n B (0)){n B (t)}) eq (46) 



Finally, substitution of this expression in eq. (38) yields 



k AB = j-^t- I dt((h B (0))(h B {t))) eq (47) 



eq JO 



By remembering the link between ensemble averages and the values associated with a single system, the equivalence 
between this equation and expression (|l2|) is then obviously proven. 
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V. GENERAL CASE 



The method of reactive flux used by Frenkel and co-workers in their simulations considers that the barrier separates 
the system into two states, A and B, and it is restricted to the evaluation of the transition rate between these initial 
and final states. In particular, it provides no information about the evolution of the intermediate configurations. In 
the previous section, we have proven that this scheme corresponds to a quasi-stationary case associated with high 
nucleation barriers. The theory we propose also allows us to evaluate the evolution and the time-dependent rate of 
change of the global crystallinity Q of the system, at any intermediate state between Q\ and Qi. 

Our scheme is formulated in terms of a Kramers-like equation accounting for the evolution of p{Q 1 t). From the 
knowledge of the number density p(Q,t) one could then evaluate the relevant quantities of the system. 

In the general case, one must work with the kinetic equation ( ^2|) which constitutes the central point of our approach. 
This equation contains the free energy barrier &(Q) and the phenomenological coefficient in the internal space L(Q). 
The barrier $(Q) could be evaluated from simulations analogous to the ones performed in reference [|. In section III 
we have obtained an explicit expression for the coefficient L(Q). However, it is necessary to transform eq. ( p7| ) into 
an expression more proper for simulations. 

With the help of eq. ( ]43| ) , one can relate the rate of change of the order parameter to the random contribution to 
the current in the following way 

(Q(t)) = [ 2 Q d6p{ 2 ,t] dQ = [ * r(Q,t)dQ (48) 
Jqi at JQi 

By substitution of the latter expression and eq. (|]|) into the Green-Kubo formula ( f27j) for £(7) (where 7 represents 
an arbitrary value of the degree of crystallization) we finally obtain 

/>oo 

k B L{i) = - / dt((Q(0))5N( 1: t)) eq (49) 
Jo 

In the equation above there appears the expected value (Q(0)), which could be identified with the initial speed of 
change Q(0) of the system which is being simulated. Additionally, the term 5N(-f,t) is given by 

5N{>y,t) = -(n 7 (t)) (50) 

where n 7 (t) = 0[Q(t) - 7]. 

With these considerations in mind, the definitive expression in a form suitable for simulation yields 

L(j) = -j^J dt(Q(0)n 7 (t)) eq (51) 

Note that the latter equation for the phenomenological coefficient is formally similar to expression 
sections II and IV. However, two significant differences exist. The first difference comes from the fact that in eq. (|12|), 
the only trajectories giving a non-vanishing contribution are the ones corresponding to configurations initially at the 
top of the barrier. The second consists of the replacement of the characteristic function ns{t) — 9[Q(t) — Qo] of 
eq. ( |l2| ) for n^(t). Therefore, by performing simulations analogous to the ones carried out by Frenkel and co-workers, 
we could obtain the phenomenological coefficient £(7) for any value 7 of the degree of crystallization. Knowing L(j), 
since the barrier $>(Q) could be evaluated from simulations, we could completely determine the Kramers-like eq. rt22) . 
By solving this equation we could determine the probability density p(Q,t) describing the dynamics of the system. 

VI. CONCLUSIONS 

In this paper we have introduced a systematic scheme to treat homogeneous nucleation, applicable to any type of 
global order parameters describing the process. This scheme is based on the internal degrees of freedom formalism 
inside the framework of non-equilibrium thermodynamics. In fact, we have provided an extension of this formalism 
to allow descriptions in terms of order parameters, relying on the construction of a non-equilibrium ensemble. 

Our treatment links theory with simulations, in the sense that the equation we propose is a Kramers-like equation as 
the ones commonly used at the theoretical level, but its coefficients are susceptible to be determined from simulations. 
Moreover, it constitutes a unifying proposal where the expressions of the method of reactive flux emerge in the 
quasi-stationary limit of a more general Kramers-type theory. 



12), analyzed in 
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The approach we have introduced provides a complete and fully consistent description of homogeneous nucleation 
for an order parameter such as the one used by Frenkel and co-workers in their simulations. Our scheme not only 
includes the reactive flux theory but also it allows a deeper study of the process, offering a framework for a theoretical 
interpretation of simulation results. We have outlined that by performing simulations analogous to the ones carried 
out in the quasi-stationary case, one could in principle obtain more complete information about the process. For 
example, we could study time-dependent nucleation and follow the evolution and the rate of change of intermediates 
degrees of crystallinity Q. 

Additionally, our scheme casts light on the drawbacks exposed by Frenkel et al. Our picture overcomes the first 
objection, in the sense that our equation constitutes the proof of the validity of a Kramers-like diffusive equation for 
global order parameters. The equation we propose also takes into account the possible dependences on Q and even 
on t of the diffusion coefficient. Moreover, it gives explicit expressions to evaluate it, which in fact are very similar to 
the ones used by Frenkel and co-workers in their simulations. 

Finally, with respect to their last concern we have to remark that the fact of working with a global description, where 
one characterizes the macroscopic state of the system by means of a single coordinate, implies that hydrodynamic 
effects, that are of local character, cannot be taken into account. However, by employing a local description of the 
process (like the one used in ref. Eq), the formalism we have introduced enables us to study hydrodynamic effects in 
nucleation processes. In order to analyze these effects, we must increase the number of variables to achieve a more 
detailed description of the system, and consider additional contributions in the Gibbs eq. ([Ti]). The way to proceed 
has been outlined in references |lj and [l8| and will be applied to nucleation processes in a future contribution. 
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